#code for DTU analysis can be found in docs DTU.Rmd. Only follow up analyses are showed in html.
rm(sumExp)
# summarise
# empirical FDR < 0.05
map_df(all_res, ~{group_by(.x, FDR = empirical_FDR < 0.05) %>% tally() %>% filter(FDR == TRUE) }, .id = "contrast")# regular FDR < 0.05
map_df(all_res, ~{group_by(.x, FDR = regular_FDR < 0.05) %>% tally() %>% filter(FDR == TRUE) }, .id = "contrast")dtu_genes <- map(all_res, ~{ .x %>% filter(empirical_FDR < 0.05) %>% pull(gene) })
#UpSetR::upset(fromList(dtu_genes),nintersects = NA)
UpSetR::upset(fromList(dtu_genes),nintersects = NA,nsets = 7)sig_res <- map_df(all_res, ~{.x %>% filter(empirical_FDR < 0.05, abs(estimates) > 1)}, .id = "contrast")
age_res <- filter(all_res$age, empirical_FDR < 0.05)
stimulation_res <- map_df(all_res, ~{.x %>% filter(empirical_FDR < 0.05, abs(estimates) > 1)}, .id = "contrast" ) %>% filter(contrast != "age")all_res %>%
map_df( ~{.x %>% filter(pval < 0.05)}, .id = "contrast") %>%
filter(!contrast %in% "age") %>%
ggplot(aes( x = estimates, y = -log10(pval), colour = empirical_FDR < 0.05 )) +
geom_point(size = 1) + scale_colour_manual(values = c("black", "red")) +
facet_wrap(~contrast, scales = "free") +
xlim(-4,4) +
theme_bw()## Warning: Removed 4 rows containing missing values (geom_point).
# age only
all_res %>%
map_df( ~{.x %>% filter(pval < 0.05)}, .id = "contrast") %>%
filter(contrast == "age") %>%
ggplot(aes( x = estimates, y = -log10(pval), colour = empirical_FDR < 0.5 )) +
geom_point(size = 1) + scale_colour_manual(values = c("black", "red")) +
facet_wrap(~contrast, scales = "free") +
xlim(-0.1,0.1) +
theme_bw()## Warning: Removed 6727 rows containing missing values (geom_point).
length(which(all_res$LPS$empirical_FDR < 0.05))## [1] 79
length(which(all_res$LPS$empirical_FDR < 0.10))## [1] 119
length(which(all_res$LPS$empirical_FDR < 0.15))## [1] 168
LPS_t <- subset(all_res$LPS, empirical_FDR < 0.05)
res_name_LPS_top = LPS_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_LPS_top)length(which(all_res$IFNy$empirical_FDR < 0.05))## [1] 65
length(which(all_res$IFNy$empirical_FDR < 0.10))## [1] 96
length(which(all_res$IFNy$empirical_FDR < 0.15))## [1] 111
IFNy_t <- subset(all_res$IFNy, empirical_FDR < 0.05)
res_name_IFNy_top = IFNy_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_IFNy_top)length(which(all_res$R848$empirical_FDR < 0.05))## [1] 59
length(which(all_res$R848$empirical_FDR < 0.10))## [1] 80
length(which(all_res$R848$empirical_FDR < 0.15))## [1] 115
R848_t <- subset(all_res$R848, empirical_FDR < 0.05)
res_name_R848_top = R848_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_R848_top)length(which(all_res$TNFa$empirical_FDR < 0.05))## [1] 17
length(which(all_res$TNFa$empirical_FDR < 0.10))## [1] 29
length(which(all_res$TNFa$empirical_FDR < 0.15))## [1] 83
TNFa_t <- subset(all_res$TNFa, empirical_FDR < 0.05)
res_name_TNFa_top = TNFa_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_TNFa_top)length(which(all_res$DEX$empirical_FDR < 0.05))## [1] 2
length(which(all_res$DEX$empirical_FDR < 0.10))## [1] 3
length(which(all_res$DEX$empirical_FDR < 0.15))## [1] 5
DEX_t <- subset(all_res$DEX, empirical_FDR < 0.05)
res_name_DEX_top = DEX_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_DEX_top)length(which(all_res$IL4$empirical_FDR < 0.05))## [1] 48
length(which(all_res$IL4$empirical_FDR < 0.10))## [1] 81
length(which(all_res$IL4$empirical_FDR < 0.15))## [1] 147
IL4_t <- subset(all_res$IL4, empirical_FDR < 0.05)
res_name_IL4_top = IL4_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_IL4_top)length(which(all_res$ATP$empirical_FDR < 0.05))## [1] 509
length(which(all_res$ATP$empirical_FDR < 0.10))## [1] 734
length(which(all_res$ATP$empirical_FDR < 0.15))## [1] 975
ATP_t <- subset(all_res$ATP, empirical_FDR < 0.05)
res_name_ATP_top = ATP_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_ATP_top)length(which(all_res$age$empirical_FDR < 0.05))## [1] 34
length(which(all_res$age$empirical_FDR < 0.10))## [1] 108
length(which(all_res$age$empirical_FDR < 0.15))## [1] 199
age_t <- subset(all_res$age, empirical_FDR < 0.05)
res_name_age_top = age_t[, c("isoform_id", "gene_id", "gene", "estimates", "se", "df", "pval", "regular_FDR", "empirical_pval", "empirical_FDR")]
createDT(res_name_age_top)